library(dplyr)
library(readxl)
library(tidyverse)
# nie<-read_xlsx("Data-NIE and WANG(1).xlsx",sheet = 1)
# pra<-read_xlsx("Mobile buttons.xlsx",sheet = 1)
# lbw<-read_xlsx("Data Lbw.xlsx",sheet = 1)
# alldata<-bind_rows(nie,pra,lbw)
# write.csv(alldata,row.names=F,"alldata.csv")
alldata<-read_csv("alldata.csv")
# alldata[alldata$Species=="Rhododendron leptothrium","Species"] <- "Rhododendron_leptothrium"
# write.csv(alldata,"alldata.csv",row.names=F)
summary(alldata)
## ForestType No.Plot Elevation Slope
## Length:104 Min. :1.00 Min. :500 Min. :0.500
## Class :character 1st Qu.:1.00 1st Qu.:500 1st Qu.:1.600
## Mode :character Median :2.00 Median :600 Median :4.000
## Mean :2.01 Mean :675 Mean :3.945
## 3rd Qu.:3.00 3rd Qu.:800 3rd Qu.:6.000
## Max. :3.00 Max. :900 Max. :8.200
## DistRiver CoordinateX CoordinateY Species
## Min. : 1.00 Min. : 1.000 Min. :0.0500 Length:104
## 1st Qu.: 6.80 1st Qu.: 2.250 1st Qu.:0.1000 Class :character
## Median :11.00 Median : 4.000 Median :0.2125 Mode :character
## Mean :10.35 Mean : 4.894 Mean :0.2719
## 3rd Qu.:13.00 3rd Qu.: 8.000 3rd Qu.:0.4188
## Max. :21.00 Max. :10.000 Max. :0.5250
## Number plotname
## Min. : 1.000 Length:104
## 1st Qu.: 1.000 Class :character
## Median : 2.000 Mode :character
## Mean : 3.298
## 3rd Qu.: 4.000
## Max. :16.000
alldata <- alldata %>% mutate(ForestType=as.factor(ForestType),Species=as.factor(Species),Elevation=as.factor(Elevation),plotname=as.factor(paste0(ForestType,No.Plot,Elevation)))
summary(alldata)
## ForestType No.Plot Elevation Slope DistRiver
## Primary :55 Min. :1.00 500:32 Min. :0.500 Min. : 1.00
## Secondary:49 1st Qu.:1.00 600:28 1st Qu.:1.600 1st Qu.: 6.80
## Median :2.00 800:22 Median :4.000 Median :11.00
## Mean :2.01 900:22 Mean :3.945 Mean :10.35
## 3rd Qu.:3.00 3rd Qu.:6.000 3rd Qu.:13.00
## Max. :3.00 Max. :8.200 Max. :21.00
##
## CoordinateX CoordinateY Species
## Min. : 1.000 Min. :0.0500 Rhododendron_decorum :19
## 1st Qu.: 2.250 1st Qu.:0.1000 Itoa_orientalis :15
## Median : 4.000 Median :0.2125 Rhododendron_leptothrium:13
## Mean : 4.894 Mean :0.2719 Ficus_heteromorpha :10
## 3rd Qu.: 8.000 3rd Qu.:0.4188 Illicium_macranthum : 6
## Max. :10.000 Max. :0.5250 Manglietia_insignis : 5
## (Other) :36
## Number plotname
## Min. : 1.000 Secondary1500: 8
## 1st Qu.: 1.000 Primary2600 : 6
## Median : 2.000 Primary3600 : 6
## Mean : 3.298 Primary1500 : 5
## 3rd Qu.: 4.000 Primary2900 : 5
## Max. :16.000 Primary3500 : 5
## (Other) :69
length(unique(alldata$Species))
## [1] 22
spabun<-as.data.frame(matrix(0,nrow=length(unique(alldata$plotname)),ncol=length(unique(alldata$Species))))
rownames(spabun)<-unique(alldata$plotname)
colnames(spabun)<-unique(alldata$Species)
summary(spabun)
## Rhododendron_leptothrium Rhododendron_decorum Ficus_heteromorpha
## Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0
## Osbeckia_mepalensis Vaccinium_duclouxii Itoa_orientalis Illicium_macranthum
## Min. :0 Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0 Max. :0
## Manglietia_insignis Beilschmiedia_robusta Sterculia_nobilis
## Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0
## Saurauia_napaulensis Padus_napaulensis Craibiodendron_yunnanense
## Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0
## Rhododendron_irroratum Ilex_corallina Lithocarpus_hancei Schefflera_fengii
## Min. :0 Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0 Max. :0
## Skimmia_arborescens Litsea_chunii Claoxylon_khasianum Schefflera_shweliensis
## Min. :0 Min. :0 Min. :0 Min. :0
## 1st Qu.:0 1st Qu.:0 1st Qu.:0 1st Qu.:0
## Median :0 Median :0 Median :0 Median :0
## Mean :0 Mean :0 Mean :0 Mean :0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0 3rd Qu.:0
## Max. :0 Max. :0 Max. :0 Max. :0
## Camellia_forrestii
## Min. :0
## 1st Qu.:0
## Median :0
## Mean :0
## 3rd Qu.:0
## Max. :0
for(i in 1:nrow(alldata)){
spabun[as.character(alldata$plotname[i]),as.character(alldata$Species[i])]<-alldata$Number[i]
}
# write.csv(spabun,"SpeciesAbundance.csv")
We altogether have 24 sites, and 22 species. Thus the species abundance matrix consists 24 rows and 22 columns.
library(tidyverse)
trait <- read_csv("dummy_trait.csv")
summary(trait)
## sp LMA LL Amass
## Length:77 Min. : 17.38 Min. : 5.38 Min. : 11.35
## Class :character 1st Qu.: 58.88 1st Qu.: 21.88 1st Qu.: 54.83
## Mode :character Median : 97.72 Median : 40.30 Median : 97.93
## Mean :132.88 Mean : 59.03 Mean :140.15
## 3rd Qu.:158.49 3rd Qu.: 81.85 3rd Qu.:174.55
## Max. :602.56 Max. :267.95 Max. :985.05
## Rmass Nmass Pmass WD
## Min. : 1.29 Min. : 0.300 Min. :0.0100 Min. :0.3100
## 1st Qu.: 7.23 1st Qu.: 0.960 1st Qu.:0.0500 1st Qu.:0.4400
## Median :13.94 Median : 1.530 Median :0.1100 Median :0.4900
## Mean :18.02 Mean : 2.744 Mean :0.1481 Mean :0.5018
## 3rd Qu.:23.80 3rd Qu.: 4.220 3rd Qu.:0.1600 3rd Qu.:0.5400
## Max. :67.58 Max. :22.690 Max. :0.7400 Max. :0.8600
## SM
## Min. :0.060
## 1st Qu.:0.550
## Median :1.350
## Mean :1.933
## 3rd Qu.:2.790
## Max. :7.980
Trait columns meaning: LMA stands for Leaf mass per area (g m^2), LL stands for Leaf lifespans (longevity, months), Amass stands for Maximum photosynthetic rates per unit mass (nnmol g^-1 s^-1), Rmass stands for Dark respiration rates per unit mass (nnmol g^-1 s^-1), Nmass stands for Leaf nitrogen per unit mass (%), Pmass stands for Leaf phosphorus per unit mass (%), WD stands for Wood density (g cm^-3), SM stands for Seed dry mass (mg).
trait %>%
pivot_longer(LMA:SM, names_to = "trait") %>%
ggplot() + aes(x=value) +
geom_histogram(position = "identity") +
facet_wrap(~ trait, scale = "free")
# transforamtion
trait2 <- trait |>
mutate(logLMA = log(LMA),
logLL = log(LL),
logAmass = log(Amass),
logRmass = log(Rmass),
logNmass = log(Nmass),
logPmass = log(Pmass),
logSM = log(SM)) |>
dplyr::select(sp, logLMA, logLL, logAmass, logRmass, logNmass, logPmass, WD, logSM)
trait2 %>%
pivot_longer(logLMA:logSM, names_to = "trait") %>%
ggplot() + aes(x=value) +
geom_histogram(position = "identity") +
facet_wrap(~ trait, scale = "free")
# much normal now
diversity <- data.frame(plot=rownames(spabun))
diversity$shannon<-vegan::diversity(spabun, index = "shannon")
diversity$No.ind <- apply(spabun,1,sum)
library(picante)
library(FD)
phylo <- read.tree("dummy_tree.newick")
plot(phylo)
diversity[,c("PD","richness")] <- pd(spabun, phylo)
diversity[,"MNTD"] <-mntd(spabun, cophenetic(phylo))
diversity[,"MPD"] <- mpd(spabun, cophenetic(phylo))
trait_mat0 <- as.matrix(trait2[, -1])
rownames(trait_mat0) <- trait2$sp
# scale
trait_mat <- scale(trait_mat0)
tmp <- trait2 %>%
filter(sp %in% colnames(spabun))
res_fd <- dbFD(trait_mat[colnames(spabun), ], spabun)
## FRic: Dimensionality reduction was required. The last 6 PCoA axes (out of 8 in total) were removed.
## FRic: Quality of the reduced-space representation = 0.8160306
diversity$FRic <- res_fd$FRic
diversity$FDiv <- res_fd$FDiv
diversity$FDis <- res_fd$FDis
diversity$RaoQ <- res_fd$RaoQ
diversity$FEVe <- res_fd$FEve
Now I have calculated all
library(stringr)
diversity$Elevation<-as.integer(str_sub(diversity$plot,start = str_length(diversity$plot)-2, end=-1))
diversity$forest <- as.factor(str_sub(diversity$plot,start = 1, end=str_length(diversity$plot)-4))
divcor <- cor(diversity[,2:12], method="pearson")
corrplot::corrplot(divcor)
Didn’t remove anything
diversity_long <- diversity %>% pivot_longer(shannon:RaoQ, names_to = "diversityindices") %>% mutate(diversityindices=factor(diversityindices, levels=c("No.ind","richness","shannon","PD","MPD","MNTD","FRic","FEVe","FDiv","FDis","RaoQ")))
ggplot(diversity_long, aes(x = Elevation, y = value, colour = forest)) + geom_point() + geom_smooth(method = lm) +facet_wrap(~diversityindices,scales = "free",ncol=3) + theme_classic()
library(lme4)
library(lmerTest)
library(performance)
ggplot(diversity_long,aes(x=value))+geom_histogram()+facet_wrap(~diversityindices,scale="free",ncol=3)
diversity %>% mutate(MNTD=log(MNTD),FRic=log(FRic),FEVe=(FEVe)^2,FDiv=(FDiv)^2,FDis=log(FDis),RaoQ=log(RaoQ)) %>%pivot_longer(shannon:RaoQ, names_to = "diversityindices") %>% mutate(diversityindices=factor(diversityindices, levels=c("No.ind","richness","shannon","PD","MPD","MNTD","FRic","FEVe","FDiv","FDis","RaoQ"))) %>% ggplot()+ aes(x = value) + geom_histogram() +facet_wrap(~diversityindices,scales = "free",ncol=3) + theme_classic()
#good enough now
diversity.trans<-diversity %>% mutate(MNTD=log(MNTD),FRic=log(FRic),FEVe=(FEVe)^2,FDiv=(FDiv)^2,FDis=log(FDis),RaoQ=log(RaoQ))
diversity.trans %>% pivot_longer(shannon:RaoQ, names_to = "diversityindices") %>% mutate(diversityindices=factor(diversityindices, levels=c("No.ind","richness","shannon","PD","MPD","MNTD","FRic","FEVe","FDiv","FDis","RaoQ"))) %>%ggplot()+aes(x = Elevation, y = value, colour = forest) + geom_point() + geom_smooth(method = lm) +facet_wrap(~diversityindices,scales = "free",ncol=3) + theme_classic()
diversity.trans$Elevation.sq <- (diversity.trans$Elevation)^2
##No. ind
indlm1 <- lm(No.ind~Elevation*forest, data=diversity.trans)
summary(indlm1) #signidicant
##
## Call:
## lm(formula = No.ind ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.083 -1.371 -0.725 1.933 5.133
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -13.916667 3.789259 -3.673 0.00151 **
## Elevation 0.045000 0.005280 8.522 4.33e-08 ***
## forestSecondary 10.450000 5.358822 1.950 0.06533 .
## Elevation:forestSecondary -0.024333 0.007467 -3.259 0.00393 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.892 on 20 degrees of freedom
## Multiple R-squared: 0.8562, Adjusted R-squared: 0.8346
## F-statistic: 39.68 on 3 and 20 DF, p-value: 1.308e-08
check_model(indlm1) # wrong!!!
indlm2 <- glm(No.ind~Elevation*forest, data=diversity.trans, family = "poisson")
summary(indlm2)
##
## Call:
## glm(formula = No.ind ~ Elevation * forest, family = "poisson",
## data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3143 -0.4419 -0.2295 0.3604 1.6520
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.9007468 0.3638483 2.476 0.0133 *
## Elevation 0.0026833 0.0004677 5.738 9.59e-09 ***
## forestSecondary 0.1031159 0.5675999 0.182 0.8558
## Elevation:forestSecondary -0.0007573 0.0007384 -1.026 0.3051
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 77.532 on 23 degrees of freedom
## Residual deviance: 12.014 on 20 degrees of freedom
## AIC: 125.53
##
## Number of Fisher Scoring iterations: 4
check_model(indlm2)
indlm3 <- glm(No.ind~Elevation, data=diversity.trans, family = "poisson")
summary(indlm3)
##
## Call:
## glm(formula = No.ind ~ Elevation, family = "poisson", data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.9752 -0.8117 -0.1612 0.9443 2.2533
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.9198876 0.2790661 3.296 0.00098 ***
## Elevation 0.0023857 0.0003615 6.600 4.11e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 77.532 on 23 degrees of freedom
## Residual deviance: 31.422 on 22 degrees of freedom
## AIC: 140.94
##
## Number of Fisher Scoring iterations: 4
check_model(indlm3)
indlm4 <- update(indlm2,.~.+Elevation.sq*forest)
summary(indlm4) #not siignificant
##
## Call:
## glm(formula = No.ind ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, family = "poisson", data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.2850 -0.4980 -0.1844 0.3629 1.5247
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.266e+00 2.316e+00 0.547 0.585
## Elevation 1.605e-03 6.766e-03 0.237 0.813
## forestSecondary 8.432e-01 3.670e+00 0.230 0.818
## Elevation.sq 7.571e-07 4.738e-06 0.160 0.873
## Elevation:forestSecondary -2.963e-03 1.077e-02 -0.275 0.783
## forestSecondary:Elevation.sq 1.560e-06 7.569e-06 0.206 0.837
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 77.532 on 23 degrees of freedom
## Residual deviance: 11.835 on 18 degrees of freedom
## AIC: 129.35
##
## Number of Fisher Scoring iterations: 4
AIC(indlm2,indlm3) # Although AIC lower, indlm3 diagnostic better
## df AIC
## indlm2 4 125.5309
## indlm3 2 140.9388
##richness
richnesslm1 <- lm(richness~Elevation*forest, data=diversity.trans)
summary(richnesslm1)
##
## Call:
## lm(formula = richness ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.7167 -0.4500 -0.2333 0.3750 2.6500
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.216667 1.315664 4.725 0.00013 ***
## Elevation -0.002333 0.001833 -1.273 0.21771
## forestSecondary 2.300000 1.860630 1.236 0.23073
## Elevation:forestSecondary -0.004000 0.002593 -1.543 0.13856
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.004 on 20 degrees of freedom
## Multiple R-squared: 0.4292, Adjusted R-squared: 0.3436
## F-statistic: 5.014 on 3 and 20 DF, p-value: 0.009405
#not significant
check_model(richnesslm1)# wrong!!!
richnesslm2 <- glm(richness~Elevation*forest, data=diversity.trans, family = "poisson")
summary(richnesslm2)
##
## Call:
## glm(formula = richness ~ Elevation * forest, family = "poisson",
## data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.8126 -0.2116 -0.1141 0.2024 1.0306
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.8761749 0.6029419 3.112 0.00186 **
## Elevation -0.0005100 0.0008551 -0.596 0.55089
## forestSecondary 0.6039792 0.8713886 0.693 0.48823
## Elevation:forestSecondary -0.0010673 0.0012608 -0.846 0.39729
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 7.7100 on 23 degrees of freedom
## Residual deviance: 4.0354 on 20 degrees of freedom
## AIC: 91.456
##
## Number of Fisher Scoring iterations: 4
#not significant
richnesslm3 <-update(richnesslm2,.~.+Elevation.sq*forest)
summary(richnesslm3)#not significant
##
## Call:
## glm(formula = richness ~ Elevation + forest + Elevation.sq +
## Elevation:forest + forest:Elevation.sq, family = "poisson",
## data = diversity.trans)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.52760 -0.30006 0.00513 0.23344 0.77838
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.298e-01 4.213e+00 0.031 0.975
## Elevation 4.763e-03 1.261e-02 0.378 0.706
## forestSecondary 6.235e+00 6.133e+00 1.017 0.309
## Elevation.sq -3.779e-06 9.014e-06 -0.419 0.675
## Elevation:forestSecondary -1.818e-02 1.847e-02 -0.984 0.325
## forestSecondary:Elevation.sq 1.231e-05 1.325e-05 0.930 0.353
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 7.7100 on 23 degrees of freedom
## Residual deviance: 3.0798 on 18 degrees of freedom
## AIC: 94.5
##
## Number of Fisher Scoring iterations: 4
##shannon
shannonlm1 <- lm(shannon~Elevation*forest,data=diversity.trans)
summary(shannonlm1)#not significant
##
## Call:
## lm(formula = shannon ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.45016 -0.15297 -0.01576 0.19303 0.43486
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.7561507 0.3216850 5.459 2.41e-05 ***
## Elevation -0.0006424 0.0004483 -1.433 0.1673
## forestSecondary 0.7578343 0.4549313 1.666 0.1113
## Elevation:forestSecondary -0.0012969 0.0006339 -2.046 0.0542 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2455 on 20 degrees of freedom
## Multiple R-squared: 0.535, Adjusted R-squared: 0.4652
## F-statistic: 7.67 on 3 and 20 DF, p-value: 0.001328
shannonlm2 <-update(shannonlm1,.~.+Elevation.sq*forest)
summary(shannonlm2)#significant
##
## Call:
## lm(formula = shannon ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.2642 -0.1451 -0.0284 0.1365 0.3963
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.046e-01 1.824e+00 -0.222 0.82697
## Elevation 5.863e-03 5.447e-03 1.076 0.29594
## forestSecondary 8.684e+00 2.580e+00 3.366 0.00344 **
## Elevation.sq -4.647e-06 3.882e-06 -1.197 0.24678
## Elevation:forestSecondary -2.516e-02 7.703e-03 -3.266 0.00429 **
## forestSecondary:Elevation.sq 1.705e-05 5.489e-06 3.105 0.00611 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2017 on 18 degrees of freedom
## Multiple R-squared: 0.7176, Adjusted R-squared: 0.6391
## F-statistic: 9.147 on 5 and 18 DF, p-value: 0.0001817
check_model(shannonlm2) #fine
##PD
PDlm1 <- lm(PD~Elevation*forest, data=diversity.trans)
summary(PDlm1)
##
## Call:
## lm(formula = PD ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.03513 -0.27121 -0.08996 0.32566 1.33305
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.5841856 0.8017674 5.718 1.35e-05 ***
## Elevation -0.0029754 0.0011172 -2.663 0.0149 *
## forestSecondary 0.5429924 1.1338703 0.479 0.6372
## Elevation:forestSecondary -0.0007405 0.0015800 -0.469 0.6444
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.6119 on 20 degrees of freedom
## Multiple R-squared: 0.476, Adjusted R-squared: 0.3973
## F-statistic: 6.055 on 3 and 20 DF, p-value: 0.004179
PDlm2 <- lm(PD~Elevation, data=diversity.trans)
summary(PDlm2)
##
## Call:
## lm(formula = PD ~ Elevation, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.05985 -0.30268 -0.09394 0.36579 1.41941
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.8556818 0.5436424 8.932 9.04e-09 ***
## Elevation -0.0033456 0.0007575 -4.416 0.000218 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.5868 on 22 degrees of freedom
## Multiple R-squared: 0.4699, Adjusted R-squared: 0.4458
## F-statistic: 19.5 on 1 and 22 DF, p-value: 0.0002182
PDlm3 <- update(PDlm1,.~.+Elevation.sq*forest)
summary(PDlm3) #not significant
##
## Call:
## lm(formula = PD ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.0358 -0.3317 -0.0536 0.3694 1.2540
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.182e+00 5.772e+00 1.244 0.229
## Elevation -1.080e-02 1.723e-02 -0.627 0.539
## forestSecondary 3.962e-01 8.162e+00 0.049 0.962
## Elevation.sq 5.587e-06 1.228e-05 0.455 0.655
## Elevation:forestSecondary -2.986e-04 2.437e-02 -0.012 0.990
## forestSecondary:Elevation.sq -3.157e-07 1.737e-05 -0.018 0.986
##
## Residual standard error: 0.6381 on 18 degrees of freedom
## Multiple R-squared: 0.4871, Adjusted R-squared: 0.3446
## F-statistic: 3.419 on 5 and 18 DF, p-value: 0.024
AIC(PDlm1,PDlm2)
## df AIC
## PDlm1 5 50.15914
## PDlm2 3 46.43283
anova(PDlm1,PDlm2)
## Analysis of Variance Table
##
## Model 1: PD ~ Elevation * forest
## Model 2: PD ~ Elevation
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 20 7.4893
## 2 22 7.5752 -2 -0.085895 0.1147 0.8922
#just keep elevation
##MPD
MPDlm1 <- lm(MPD~Elevation*forest,data=diversity.trans)
summary(MPDlm1)#not significant
##
## Call:
## lm(formula = MPD ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.95745 -0.14704 0.02182 0.11525 0.54686
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.9718876 0.4656388 4.235 0.000406 ***
## Elevation -0.0012258 0.0006489 -1.889 0.073464 .
## forestSecondary -0.1566694 0.6585127 -0.238 0.814368
## Elevation:forestSecondary 0.0003193 0.0009176 0.348 0.731529
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3554 on 20 degrees of freedom
## Multiple R-squared: 0.2228, Adjusted R-squared: 0.1062
## F-statistic: 1.911 on 3 and 20 DF, p-value: 0.1604
MPDlm2 <- update(MPDlm1,.~.+Elevation.sq*forest)
summary(MPDlm2)#not significant
##
## Call:
## lm(formula = MPD ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.03991 -0.10106 0.01737 0.12620 0.52058
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.601e+00 3.274e+00 1.406 0.177
## Elevation -9.142e-03 9.774e-03 -0.935 0.362
## forestSecondary -5.342e+00 4.630e+00 -1.154 0.264
## Elevation.sq 5.654e-06 6.966e-06 0.812 0.428
## Elevation:forestSecondary 1.593e-02 1.382e-02 1.153 0.264
## forestSecondary:Elevation.sq -1.115e-05 9.851e-06 -1.132 0.272
##
## Residual standard error: 0.3619 on 18 degrees of freedom
## Multiple R-squared: 0.2744, Adjusted R-squared: 0.0729
## F-statistic: 1.362 on 5 and 18 DF, p-value: 0.2846
##MNTD
MNTDlm1 <- lm(MNTD~Elevation*forest,data=diversity.trans)
summary(MNTDlm1)#not significant
##
## Call:
## lm(formula = MNTD ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.70494 -0.15764 0.04962 0.19746 0.91350
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.6952460 0.7429117 0.936 0.3605
## Elevation -0.0018248 0.0010352 -1.763 0.0932 .
## forestSecondary -0.2976132 1.0506358 -0.283 0.7799
## Elevation:forestSecondary 0.0006542 0.0014640 0.447 0.6598
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.567 on 20 degrees of freedom
## Multiple R-squared: 0.1957, Adjusted R-squared: 0.07502
## F-statistic: 1.622 on 3 and 20 DF, p-value: 0.2159
MNTDlm2 <-update(MNTDlm1,.~.+Elevation.sq*forest)
summary(MNTDlm2)#not significant
##
## Call:
## lm(formula = MNTD ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.81234 -0.12855 0.02618 0.20882 0.80610
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.823e+00 5.206e+00 1.118 0.278
## Elevation -1.726e-02 1.554e-02 -1.111 0.281
## forestSecondary -8.755e+00 7.362e+00 -1.189 0.250
## Elevation.sq 1.103e-05 1.108e-05 0.995 0.333
## Elevation:forestSecondary 2.612e-02 2.198e-02 1.188 0.250
## forestSecondary:Elevation.sq -1.819e-05 1.567e-05 -1.161 0.261
##
## Residual standard error: 0.5756 on 18 degrees of freedom
## Multiple R-squared: 0.254, Adjusted R-squared: 0.04684
## F-statistic: 1.226 on 5 and 18 DF, p-value: 0.3374
##FRic
FRiclm1 <- lm(FRic~Elevation*forest,data=diversity.trans)
summary(FRiclm1) #not significant
##
## Call:
## lm(formula = FRic ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.8991 -0.1224 0.1448 0.5541 1.0231
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.5044259 1.1066559 1.359 0.189
## Elevation -0.0002271 0.0015421 -0.147 0.884
## forestSecondary 0.9192560 1.5650477 0.587 0.564
## Elevation:forestSecondary -0.0011251 0.0021808 -0.516 0.612
##
## Residual standard error: 0.8446 on 20 degrees of freedom
## Multiple R-squared: 0.04473, Adjusted R-squared: -0.09856
## F-statistic: 0.3121 on 3 and 20 DF, p-value: 0.8164
FRiclm2 <- update(FRiclm1, .~.+Elevation.sq*forest)
summary(FRiclm2) # significant
##
## Call:
## lm(formula = FRic ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.4962 -0.3510 0.1637 0.4185 1.0458
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.099e+01 7.031e+00 -1.563 0.1355
## Elevation 3.738e-02 2.099e-02 1.781 0.0918 .
## forestSecondary 2.416e+01 9.943e+00 2.430 0.0258 *
## Elevation.sq -2.686e-05 1.496e-05 -1.796 0.0893 .
## Elevation:forestSecondary -7.110e-02 2.969e-02 -2.395 0.0277 *
## forestSecondary:Elevation.sq 4.998e-05 2.116e-05 2.362 0.0296 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7774 on 18 degrees of freedom
## Multiple R-squared: 0.2718, Adjusted R-squared: 0.06949
## F-statistic: 1.344 on 5 and 18 DF, p-value: 0.2912
check_model(FRiclm2) # diagnostic not good
##FEVe
FEVelm1 <- lm(FEVe~Elevation*forest,data=diversity.trans)
summary(FEVelm1)
##
## Call:
## lm(formula = FEVe ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.29327 -0.09582 0.01052 0.12074 0.29076
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.7660142 0.2237837 3.423 0.00269 **
## Elevation -0.0003705 0.0003118 -1.188 0.24868
## forestSecondary 0.2298707 0.3164779 0.726 0.47605
## Elevation:forestSecondary -0.0004105 0.0004410 -0.931 0.36299
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1708 on 20 degrees of freedom
## Multiple R-squared: 0.2949, Adjusted R-squared: 0.1892
## F-statistic: 2.789 on 3 and 20 DF, p-value: 0.06712
#not significant
FEVelm2<-update(FEVelm1, .~.+Elevation.sq*forest)
summary(FEVelm2)# not significant
##
## Call:
## lm(formula = FEVe ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.253418 -0.119571 0.005495 0.110520 0.283927
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.693e-01 1.601e+00 -0.293 0.773
## Elevation 3.349e-03 4.779e-03 0.701 0.492
## forestSecondary 1.254e+00 2.264e+00 0.554 0.587
## Elevation.sq -2.657e-06 3.406e-06 -0.780 0.446
## Elevation:forestSecondary -3.493e-03 6.759e-03 -0.517 0.612
## forestSecondary:Elevation.sq 2.201e-06 4.817e-06 0.457 0.653
##
## Residual standard error: 0.177 on 18 degrees of freedom
## Multiple R-squared: 0.3186, Adjusted R-squared: 0.1294
## F-statistic: 1.683 on 5 and 18 DF, p-value: 0.1896
##FDiv
FDivlm1 <- lm(FDiv~Elevation*forest,data=diversity.trans)
summary(FDivlm1)# significant
##
## Call:
## lm(formula = FDiv ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.30058 -0.08177 0.04738 0.10321 0.15507
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.4506895 0.1770097 2.546 0.0192 *
## Elevation 0.0002208 0.0002467 0.895 0.3814
## forestSecondary 0.5176315 0.2503295 2.068 0.0518 .
## Elevation:forestSecondary -0.0007461 0.0003488 -2.139 0.0450 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1351 on 20 degrees of freedom
## Multiple R-squared: 0.2109, Adjusted R-squared: 0.09249
## F-statistic: 1.781 on 3 and 20 DF, p-value: 0.1831
check_model(FDivlm1) # not good enough but make sense
FDivlm2 <- update(FDivlm1,.~.+Elevation.sq*forest)
summary(FDivlm2) # not significant
##
## Call:
## lm(formula = FDiv ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.25571 -0.07467 0.01305 0.09138 0.18628
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -9.403e-01 1.235e+00 -0.761 0.456
## Elevation 4.409e-03 3.688e-03 1.195 0.247
## forestSecondary 1.265e+00 1.747e+00 0.724 0.478
## Elevation.sq -2.991e-06 2.628e-06 -1.138 0.270
## Elevation:forestSecondary -2.996e-03 5.216e-03 -0.574 0.573
## forestSecondary:Elevation.sq 1.607e-06 3.717e-06 0.432 0.671
##
## Residual standard error: 0.1366 on 18 degrees of freedom
## Multiple R-squared: 0.2743, Adjusted R-squared: 0.07268
## F-statistic: 1.361 on 5 and 18 DF, p-value: 0.285
##FDis
FDislm1 <- lm(FDis~Elevation*forest,data=diversity.trans)
summary(FDislm1)# significanyt
##
## Call:
## lm(formula = FDis ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.57008 -0.17965 0.04284 0.18087 0.30536
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.7423284 0.3379152 2.197 0.0400 *
## Elevation -0.0003485 0.0004709 -0.740 0.4679
## forestSecondary 1.0471963 0.4778842 2.191 0.0404 *
## Elevation:forestSecondary -0.0014660 0.0006659 -2.202 0.0396 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.2579 on 20 degrees of freedom
## Multiple R-squared: 0.4356, Adjusted R-squared: 0.351
## F-statistic: 5.146 on 3 and 20 DF, p-value: 0.008458
FDislm2 <-update(FDislm1, .~.+Elevation.sq*forest)
summary(FDislm2) # very very significant
##
## Call:
## lm(formula = FDis ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.39544 -0.09303 0.03150 0.09719 0.32101
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.229e+00 1.726e+00 -2.450 0.024762 *
## Elevation 1.462e-02 5.154e-03 2.836 0.010954 *
## forestSecondary 1.143e+01 2.441e+00 4.683 0.000185 ***
## Elevation.sq -1.069e-05 3.673e-06 -2.910 0.009334 **
## Elevation:forestSecondary -3.273e-02 7.289e-03 -4.490 0.000283 ***
## forestSecondary:Elevation.sq 2.233e-05 5.195e-06 4.299 0.000432 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.1909 on 18 degrees of freedom
## Multiple R-squared: 0.7218, Adjusted R-squared: 0.6445
## F-statistic: 9.34 on 5 and 18 DF, p-value: 0.0001599
check_model(FDislm2) #good
##RaoQ
RaoQlm1 <-lm(RaoQ~Elevation*forest,data=diversity.trans)
summary(RaoQlm1) # not significance
##
## Call:
## lm(formula = RaoQ ~ Elevation * forest, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.1537 -0.2838 0.1330 0.3700 0.5420
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1.8620844 0.6464953 2.880 0.00925 **
## Elevation -0.0009568 0.0009009 -1.062 0.30086
## forestSecondary 1.2459579 0.9142824 1.363 0.18810
## Elevation:forestSecondary -0.0015943 0.0012740 -1.251 0.22523
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4934 on 20 degrees of freedom
## Multiple R-squared: 0.3235, Adjusted R-squared: 0.222
## F-statistic: 3.188 on 3 and 20 DF, p-value: 0.04599
RaoQlm2 <- update(RaoQlm1,.~.+Elevation.sq*forest)
summary(RaoQlm2) # significant
##
## Call:
## lm(formula = RaoQ ~ Elevation + forest + Elevation.sq + Elevation:forest +
## forest:Elevation.sq, data = diversity.trans)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.85525 -0.20003 0.08947 0.22197 0.52818
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.822e+00 3.605e+00 -1.892 0.07462 .
## Elevation 2.519e-02 1.076e-02 2.340 0.03099 *
## forestSecondary 1.918e+01 5.098e+00 3.763 0.00142 **
## Elevation.sq -1.867e-05 7.670e-06 -2.435 0.02553 *
## Elevation:forestSecondary -5.559e-02 1.522e-02 -3.652 0.00182 **
## forestSecondary:Elevation.sq 3.857e-05 1.085e-05 3.556 0.00226 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3985 on 18 degrees of freedom
## Multiple R-squared: 0.6028, Adjusted R-squared: 0.4925
## F-statistic: 5.463 on 5 and 18 DF, p-value: 0.003131
check_model(RaoQlm2)# not good
Significant models: indlm3, shannonlm2, PDlm2, FDivlm1, FDislm2